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Abstract 

A singularity on the negative fugacity axis of the hard-square lattice gas is investigated in 
terms of numerical diagonalization of transfer matrices. The location of the singular point z~ 
and the critical exponent v are accurately determined by the phenomenological renormalization 
technique as —0.11933888188(1) and 0.416667(1), respectively. It is also found that the central 
charge c and the dominant scaling dimension x a are —4.399996(8) and —0.3999996(7), respec- 
tively. These results strongly support that this singularity belongs to the same universality class 
as the Yang-Lee edge singularity (c = — ^, x a — — | and v = 

1 Introduction 

For systems of hard particles, the radius of convergence of the fugacity series expansion is known to 
be generally determined by a non-physical singularity on the negative real fugacity z axis [l], ||, |3|. 
Here, let us consider the hard-square lattice gas (hard squares), which is one of the simplest models 
describing the fluid-solid transition in two dimensions. This model exhibits a second-order transition 
at 

zt = 3.796255174(3) (1) 

from the disordered phase to the \/2 x v2 phase || (uncertainty in the last decimal digits is given 
by the figure in parentheses). It is also precisely confirmed that the transition belongs to the same 
universality class as the Ising ferromagnet ||. 

On the other hand, the alternating behavior in the sign of the coefficients in the fugacity series 
of the reduced pressure 0: 
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manifests that a dominant singularity exists on the negative fugacity axis. By series analysis @, 
the location of the singularity is estimated as 

z~ = -0.1193388809(10). (3) 
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Note that the absolute magnitude of z~ is about 30 times smaller than that of z£ . This fact makes 
it quite difficult to derive precise information around the physical critical point from the fugacity 
series Q. Similar situations are commonly observed in other hard-core systems. 

This non-physical singularity is referred to as the hard-core singularity or the repulsive- core 
singularity. At this point thermodynamic quantities are known to exhibit interesting non-trivial 
behavior. For example, the reduced pressure or equivalently the reduced free-energy density behaves 
as 

- P(*) = /(*) ~ (* - (4) 
for z — > z~+. The singularity is characterized by the exponent <j). For the hard-square lattice gas 
(f> is found to be about 0.83337 by series analysis ||. 

In 1984, Poland |J investigated the fugacity series for a variety of hard-core lattice gases and 
hard-core gases in continuous space in two, three and higher dimensions, and proposed that the 
exponent <j) is universal, that is, it depends only on the dimensionality of space. Another study on 
further models by Baram and Luban jji]] also supports this conjecture. 

In addition, Lai and Fisher [11| pointed out that the hard-core singularity can be identified with 



the Yang-Lee edge singularity pa, pi 14, 15, 16, [L7| (the exponent <j> of the hard-core singularity 



relates to Fisher's exponent a of the Yang-Lee edge singularity by cj> = a + 1). It means, if it is 
true, that the singular point z = z~ can be regarded as a conventional critical point; the correlation 
length diverges and consequently ordinary finite-size scaling analysis does work. 

As far as we know, all the studies with respect to the hard-core singularity are done by means 



of series analysis up to the present time (except for the analytic solution for hard hexagons [10, |18[ 
and the trivial cases in dimensions less than two). In this paper, we investigate this singularity of 
the hard-square lattice gas by another approach — numerical diagonalization of transfer matrices 
and phenomenological renormalization analysis. 

The transfer-matrix method is one of the most powerful numerical methods in investigating 
statistical mechanical models in low dimensions [jl{J. It has some advantages such as the absence 
of statistical errors and critical slowing down, which are observed in Monte Carlo simulations. 
Moreover, in two dimensions, the conformal invariance associated with a critical system yields a 



great deal of interesting consequences [gy, 21 , |2^] . Especially, useful asymptotically exact relations 
between universal quantities, such as the critical exponent, and eigenvalues of the transfer matrix 
are given. This enables us to perform precise analysis about the criticality even in the non-unitary 
case with a negative central charge, such as the present hard-core singularity, as shown later. 

The organization of the present paper is as follows: In §2 we describe the transfer-matrix method 
and give the explicit notation of the transfer matrix for the hard-square lattice gas. In §3 the critical 
point of the hard-core singularity and some universal quantities are investigated by means of the 
phenomenological renormalization technique. In addition, in §4 we consider the critical eigenvalue 
spectrum of the transfer matrix precisely, which relates to the operator content of the corresponding 
conformal theory. We give a summary and some discussion in the final section. 

2 Transfer Matrix for Hard-Square Lattice Gas 

The partition function of the hard-square lattice gas is defined as follows: 

Z = J2^U(l-^jl (5) 

{•<} {id) 

where Sj is a one-bit binary number, which describes whether the ith site on the square lattice 
is occupied (sj = 1) or vacant (sj = 0). The product in eq. @ runs over all nearest-neighbor 
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Figure 1: Two neighboring unit slices of the row-to-row (a) and the diagonal-to-diagonal (b) transfer 
matrices. The transfer direction is chosen in the horizontal direction (right to left) in both cases. 
In the vertical direction, periodic boundary conditions are applied. 



pairs of the lattice sites. The size of hard squares is so that any two of them can not occupy on 
nearest-neighbor positions simultaneously. 

The partition function of the system ofLxM sites with periodic boundary conditions can be 
expressed in terms of the transfer matrix fl9fl : 

Z LM = Tr (T L ) M . (6) 

The transfer matrix Tl is square and of order 2 L . Note there remains some options in choosing 
the transfer direction and dividing the total Hamiltonian into a sum of slice Hamiltonians |ll| . In 
this paper, we consider the following two types of transfer matrices, which have distinct transfer 
directions with each other. 



2.1 Row-to-row transfer matrix 

The first one is the row-to-row transfer matrix. The transfer direction is chosen to be parallel to a 
set of lattice edges (Fig. [j] (a)). The (r\a) element is explicitly written as 

L L 

T t , L {r\a) = z^i=i Sl Y[(l - s iSi+1 )(l - Si ti){l - UU+i), (7) 

8=1 

where a = (si, S2, ■ ■ ■ , sl) and r = (ti, t2, ■ ■ ■ ,t£) respectively denote one of 2 L configurations on 
the two neighboring unit slices respectively (sj, ti = 1 or 0). Periodic boundary conditions along 
the unit slices (the vertical direction in Fig. [I] (a)) are applied (si + L = Sj and ij+L = ti). 

It should be noticed that quite a number of matrix elements should vanish owing to the presence 
of the hard-core interaction. Especially, the number of rows or columns containing at least one non- 
zero element is only + \/5)) = (1-61803 • • -) L for L>1. This is much smaller than 2 L . One 
can construct an algorithm for multiplication of the matrix to a trial vector almost in this highly 
restricted space as shown by Todo and Suzuki ||. Thus, dominant eigenvalues of the transfer 
matrix can be calculated effectively up to the system width L = 28, which is quite larger than that 



achieved by the standard sparse- matrix factorization |5L 23]. 
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2.2 Diagonal-to-diagonal transfer matrix 

We also consider the diagonal-to-diagonal transfer matrix. In this case, the transfer direction is 
rotated by the amount of ? to a set of lattice edges (Fig. [l] (b)). Periodic boundary conditions are 
also applied in the orthogonal direction. Matrix elements are written as follows: 

L L 

Td,U T W) = z^^ Si - - s i+1 U), (8) 

i=l 

In this case, any two particles on the sites in a unit slice do not interact with each other. Conse- 
quently, all the rows or columns of the matrix contains at least one none-zero element. By 
using the standard sparse-matrix factorization technique [||, 23], numerical matrix multiplication 
to a vector can be carried out up to L = 20. 

2.3 Eigenvalues, free-energy density and correlation length 

In the M — > oo limit, the reduced free-energy density fi is given in terms of the largest eigenvalue 
Ai of the transfer matrix as 

h = lim ( -TT7 lo & Z LM 



M^oo V LM 

= --logAi. (9) 

The inverse correlation length in the transfer direction, which is also referred to as the gap, is 
calculated as 

C L 1 = 9L = log ^ , (10) 

where A2 is the second largest eigenvalue in absolute magnitude. 

It should be noted that for z < a number of matrix elements in eqs. (|7j) and (||) become 
negative. Consequently, there is no guarantee that the largest eigenvalue is unique. In general, it 
may appear as a member of degenerating eigenvalues or one of a complex-conjugate pair. It is in a 
sharp contrast to the case for z > 0, where the largest eigenvalue is always real and non-degenerate 
as a consequence of the Perron- Frobenius theorem p4|]. 

However, in the present case, one can easily show that there exists a certain fugacity region 
(z£ < z < 0), where the largest eigenvalue is actually real and non-degenerate, as the following: 

First of all, note that all the elements of and T& l are real. Then, it is apparent coefficients of 
their characteristic polynomials are all real. This means that imaginary eigenvalues should appear 
as complex-conjugate pairs, if they exist. On the other hand, it can be easily seen that at z = 
the largest eigenvalue is unity and all the other eigenvalues are zero. Consequently, the largest 
eigenvalues is real and non-degenerate at least between the z = and the point z = z£ , where the 
largest and the second-largest eigenvalues first cross with each other [^] . For z < z£, they form a 
complex-conjugate pair. In the following, we work only in the region z£ < z < 0. 

3 Phenomenological Renormalization 

The largest and the second-largest eigenvalues of the transfer matrix are calculated by the power 



method [26| up to L = 28 and L = 20 for the row-to-row and diagonal-to-diagonal transfer matrices, 



respectively. In the present case, other improved variations of the power method such as the 
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Table 1: Critical point z L , exponent vl, scaled gap Lgi J (z)/2'K at z* L and central charge cl for the 
row-to-row transfer matrix. Their extrapolated values to L — » oo are also listed in the bottom row. 



L 




4 








Lg L /2ir 




CL 


2 


-0 


11569699529389 





421849529 





4814837782 





417414138 


3 


-0 


11856124181461 





416582977 





4394700591 





438343788 


4 


-0 


11908596327934 





415750210 





4233109453 





430386102 


5 


-0 


11923864086390 





415659639 





4150498994 





421997003 


6 


-0 


11929304185438 





415705242 





4104150242 





416258939 


7 


-0 


11931538407218 





415775276 





4076315079 





412555119 


8 


-0 


11932569455336 





415846360 





4058502112 





410099405 


9 


-0 


11933093621984 





415912883 





4046435075 





408387975 


10 


-0 


11933381954516 





415973444 





4037860694 





407135062 


11 


-0 


11933550958958 





416027918 





4031528434 





406179637 


12 


-0 


11933655239965 





416076654 





4026705275 





405427508 


13 


-0 


11933722350045 





416120179 





4022938659 





404820630 


14 


-0 


11933767079976 





416159066 





4019936087 





404321316 


15 


-0 


11933797791846 





416193868 





4017501012 





403903960 


16 


-0 


11933819423819 





416225088 





4015496991 





403550501 


17 


-0 


11933835002255 





416253176 





4013826680 





403247800 


18 


-0 


11933846442079 





416278521 





4012418994 





402986066 


19 


-0 


11933854989180 





416301464 





4011220978 





402757844 


20 


-0 


11933861474296 





416322296 





4010192489 





402557360 


21 


-0 


11933866463558 





416341268 





4009302622 





402380074 


22 


-0 


11933870350381 





416358598 





4008527257 





402222367 


23 


-0 


11933873413008 





416374471 





4007847341 





402081322 


24 


-0 


11933875851387 





416389050 





4007247649 





401954561 


25 


-0 


11933877811309 





416402474 





4006715898 





401840128 


26 


-0 


11933879400494 





416414865 





4006242084 





401736403 


27 


-0 


11933880699512 





416426329 





4005817991 





401642027 


oo 


-0.11933888189(2) 


0.416668(2) 


0.3999996(7) 


0.400000(2) 



conjugate-gradient method or the Lanczos method [27] are found to be less numerically stable. In 
the computation, we use quadratic-precision real numbers instead of double-precision numbers, so 
that finite-size quantities [z* L in eq. (|TT|) , etc.) can be obtained with the uncertainty of 10~ 20 or 
less. As seen below, finite-size data exhibit quite rapid convergence to the thermodynamic limit. 
Then, the quite high accuracy for them is required to make meaningful extrapolations. 

The critical point is estimated by solving the phenomenological renormalization equations: 

Lg L (zl)=L'g L/ (zt) (11) 

for two successive system sizes L and L' = L + 1. As shown in Table I and II, there is no alternating 
behavior in finite-size corrections with respect to the parity of L. It is in contrast with the case for 
the physical critical point at z = z^~ , where the systems with even and odd width have different 
finite-size corrections with each other, reflecting the anti-ferromagnetic-like y/2 x v2 ordering in 
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Table 2: Critical point z L , exponent vl, scaled gap LgL(z)/2ir at z* L and central charge cl for the 
diagonal-to-diagonal transfer matrix. Their extrapolated values to L — > oo are also listed in the 
bottom row. 



J 

L 




Z L 








Lg L /2ir 




CL 


2 


-0 


11876248705617 





447460528 





4308877716 





320446273 


3 


-0 


11923273442241 





431550594 


(J 


A 1 -"l / ' "1 P /I 00*7 

4126164837 





355461268 


4 


-U 


11930719600997 





425758018 


(J 


4068301638 





370796613 


5 


-0 


11932661625882 





422898114 


(J 


4042554292 





379106387 


6 


-0 


11933327955420 





421248828 





4028891193 





on <i 01 1 TO 

384181172 


I 


-0 


H933601031532 





A OAAAAArO 

420200952 





A C\C~\C\'~? C\V f\'1C\ 

4020795070 





OOTfO A rp^ 

387534567 


8 


-0 


11933727916302 





419489095 





4015617111 





389879246 


9 


-0 


11933792672807 





/II Oflfll AfO 

418981068 





>i ai rii 1 o r» o 1 

4012113921 





391590106 


1U 


-U 




n 
U 


/i i Q^n/i /lev 
41oOU44o ( 


n 
U 


4UUyooyUUo 


u 


oyzoouyou 


1 1 
1 1 


-0 


11933848910493 





418316797 





4007829378 





393881552 


12 


-0 


11933861540512 





418091545 





4006468626 





394674591 


13 


-0 


11933869564988 





417911541 





4005421319 





395314964 


14 


-0 


11933874838811 





417765191 





4004599259 





395840334 


15 


-0 


11933878407038 





417644431 





4003943056 





396277293 


16 


-0 


11933880882938 





417543502 





4003411561 





396645072 


17 


-0 


11933882639297 





417458199 





4002975557 





396957882 


18 


-0 


11933883909801 





417385385 





4002613844 





397226416 


19 


-0 


11933884844969 





417322683 





4002310753 





397458855 


oo 


-0.11933888188(1) 


0.416667(1) 


0.3999994(9) 


0.399999(1) 



the ordered phase ||. Then, iterated fits || of 3rd order yields 

z~ = -0.11933888189(2) (12) 
for the row-to-row transfer matrix and 

z - = -0.11933888188(1) (13) 

for the diagonal-to-diagonal one. They agree excellently with each other. They are also consistent 
with eq. but much accurate by about two orders of magnitude. 
The exponent v can be found from derivatives of the gap: 

at z = z* L . The results are also listed in Table I and II. They are extrapolated to the thermodynamic 
limit by iterated fits as 

v = 0.416668(2) (15) 

and 

v = 0.416667(1) (16) 

in the row-to-row and diagonal-to-diagonal cases, respectively. These results indicate v = J^, which 
is expected to be exact for the Yang-Lee edge singularity in two dimensions [16]. 
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At the critical point, the scaled gap Lgi,{z) in eq. (11) is believed to be universal [28, 2S], that 



is, 



2irx r 



(17) 



for L — * oo, where x a would be the scaling dimension of the dominant scaling operator and relate 
to the exponent of the correlation function r] by r\ = 1x a . It is also expected that the coefficient of 



the leading finite-size correction to the critical free-energy density becomes universal 3C, 31]: 

7TC 



fo 



6L 2 



+ 



(18) 



where fo is the critical free-energy density of the bulk, and c is the central charge of the Virasoro 



algebra, which characterize the critical theory of the model [20]. 
In practice, finite-size estimates for x a and c are calculated by 

1 



Xa,L = —Lg L (z* L ) 



and 



6 L L 



2 r'2 



cl 



7T (L' 2 - L 2 ) 



{fu(zl)-h(zl)}, 



(19) 



(20) 



respectively. As seen in Table || and ||, both of x^l and cl tend to rapidly converge to the same 

c = 0.399999(1) (21) 



value | for L — > oo. The best estimates are 



and 



x (j = 0.3999996(7), 



(22) 



respectively. Following these results, one may conclude that the singularity at z = z~ would 
be described by a conformal field theory with c = | and x a = |. However, this theory does 
not appear in the classification of the unitary minimal series by Friedan, Qiu and Shenker [21] 
(c = 1 — 6/m(m + 1) with m = 3, 4, 5, • ■ •). Is it an absolutely new class of universality? 

Another interpretation for our numerical results are possible. Suppose the largest eigenvalue Ai 
corresponds not to the identity operator but to the primary operator with a negative dimension. 
Then, the second largest one A2 now corresponds to the identity operator with the dimension 0. 
Accepting this assumption, one should calculate the true central charge and scaling dimension (say 
c and x a , respectively) same as eqs. (pi) and (M), but using the modified free-energy density and 
gap @: 



and 



9L = log 



log A 2 

A2 
Ai 



(23) 



(24) 



instead of those defined in eqs. (S) and (|l0|). In the expressions the largest and second-largest 
eigenvalues, Ai and A2 change places with each other in comparison with the former definitions. 
Thus, c and related with c and x a above mentioned as 



c — 12xo- 
-4.399996(8) 



(25) 
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Table 3: Eigenvalue spectrum of the transfer matrix at z = z~ . The degeneracy of each level is 
denoted by d in the second column. The degeneracy and the scaling dimensions of the operators 
appearing in the modular invariant partition function for the Yang-Lee edge singularity ]l7J are 
also listed in the fourth and fifth columns, respectively. 







Present result 




Reference [|17j] 


k 


d 


X'z- 


d 


X u 


1 


1 


-0.3999996(7) 


1 


•2 

5 


2 


1 





1 





3 4 


2 


0.6001(2) 


2 


-1 + 1 


5 


1 


1.599(2) 


3 


-| + 2 


6 7 


2 


1.598(3) 






8 9 


2 


2.0001(1) 


2 


+ 2 


10 11 


2 


2.607(8) 


4 


-| + 3 


12 13 


2 


2.606(9) 






14 15 


2 


2.996(8) 


2 


+ 3 



and 

Xq- X(j 

= -0.3999996(7), (26) 

respectively. These values coincides quite excellently with Cardy's conjecture that the Yang-Lee 
edge singularity would be described by the non-unitary conformal field theory with c = — ^ and 
x a = — |. It is also confirmed that v = ^ and r\ = 2x a = — | satisfy Fisher's scaling relation |l5|] 
i/ = 2/(c£ + 2 — 77), where d is the dimensionality of space. In the next section, we investigate the 
critical eigenvalue spectrum of the transfer matrix in order to confirm our tentative conclusion more 
precisely. 



4 Eigenvalue Spectrum of Transfer Matrix at Critical Point 



The complete spectrum of the transfer matrix at the critical point contains more information 
about the criticality. It is known that at the critical point, there is asymptotically one-to-one 
correspondence between the eigenvalues of the transfer matrix and the operators appearing in the 
corresponding conformal field theory |^, 29]. 

We calculated all the eigenvalues of the row-to-row transfer matrix (R) at z = z~ by the 
combination of the Householder transformation and the QR algorithm J26|. The corresponding 
scaling dimensions are then calculated as 



S k = lim ±.L§Pte 



with 



9 L = lo S 



A2 
A*. 



(27) 
(28) 
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where is the fcth-largest eigenvalue. Note that the gap in eq. fl28| ) is measured from the second- 
largest eigenvalue A2 following the assumption in the last section. 

The extrapolated values for the first 15 eigenvalues (k = 1, 2, • • • , 15) are listed in Table III with 
their degeneracy. (The results for x\ has been obtained as x a in the last section (eq. (|26|)), and x^ 
is apparently by definition (eq. (|28|)).) 

First, it should be noted that the dimension 1 is absent in the eigenvalue spectrum (Table III). 
This supports our conclusion, because the level-1 descendants of the identity operator generally 
vanish [20|. On the other hand, the dimension — | + 1 is present, which corresponds to the level-1 
secondary operators of the operator with the dimension x\ = x a . 

Moreover, the whole eigenvalue spectrum, at least up to k = 15, coincides completely with the 
operator content of the modular invariant partition function for c = — ^ conformal theory obtained 
by Itzykson et al. |jT7j| , including the degeneracy. (Their results are also listed in the fourth and 
fifth columns of Table III). 



5 Summary and Discussion 

In the present paper, we investigated the hard-core singularity of the hard-square lattice gas by the 
transfer-matrix method and the phenomenological renormalization. In special, for the row-to-row 
transfer matrix we perform diagonalization up to the quite large system width L = 28 by using the 
improved sparse-matrix factorization technique. Finite-size analysis for the two different transfer 
matrices yields the highly consistent results with each other. In terms of the phenomenological 
renormalization technique, we obtained the accurate estimate for the critical point z c , which is 
consistent with the previous studies but much accurate. 

The exponent the central charge c and the dominant scaling dimension x a are also determined 
accurately. The results indicate strongly that the hard-core singularity is described by the non- 
unitary conformal field theory with c = — in other words, it belongs to the same universality 
class as the Yang-Lee edge singularity. 

As far as we know, the present result for the critical point is one of the most accurate numerical 
estimates by means of the computer simulation for the non-trivial problems in statistical physics. 
The rapidness of convergence or the smallness of finite-size corrections would be due to the extreme 
simplicity of the critical theory — there are only two primary operators (x = and — |). 

Finally, we comment that the hard-square lattice gas can be considered as a special limit of the 
Ising anti-ferromagnet on the square lattice fl3| {(3 J — > —00 and j3H — > 00 with [3{H + 2 J) fixed, 
where (3 is the inverse temperature, J the nearest-neighbor anti-ferromagnetic interaction (J < 0) 
and the H is the magnetic field). In the latter model, the Yang-Lee edges are expected to form 
a couple of critical lines in the complex H plane above the critical temperature. The connection 
between them is an interesting problem, but it remains as an open question. 
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